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ABSTRACT 

Using two dimensional simulations, we study the formation of structures with a high- 
density contrast by magnetohydrodynamic waves in regions in which the ratio of ther- 
mal to magnetic pressure is small. The initial state is a uniform background perturbed 
by fast-mode wave. Our most significant result is that dense structures persist for far 
longer in a two-dimensional simulation than in the one-dimensional case. Once formed, 
these structures persist as long as the fast-mode amplitude remains high. 
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1 INTRODUCTION 

One of the outstanding problems in star formation is the 
determination of the mechanisms that are responsible for 
the observed dumpiness in star forming regions. Since 
the ratio, f3, of the thermal to ma gnetic pressure is ob- 
served to be small Scratched Il99glh it is natural to sup- 
pose that magnetohydrodynamic (MHD) waves play a cru- 
cial role. For this reason, there have been a number of 
three-dimensional (3D) simulations of the effect of MHD 
waves on an isothermal plasma in which is small (e.g. 
Ballcstcros-Parcdes & MacLow 2002; Padaan & Nordlunc 
20021: iGammie et al.1 120031: iLi et alJ l2004t iNakamura fc L 
2003 : T^zoti^^emad^niet^lT 1200,4) . These do, indeed, 
show that density inhomogeneities can be formed in this way 
with statistical properties that are consistent with the ob- 
servations. However, these calculations generally contain so 
many ingredients that they obscure the fundamental mech- 



In two previo us papers (iFalle fc Hartauistll 2002 , here- 
after FH02 and iLim. Falle fc HartauislJ 12005, hereafter 
LFH05), we carried out ID calculations of MHD waves in 
a low-/3 plasma that were designed to elucidate the mech- 
anism by which clumps are formed. In FH02, we showed 
that clumps with a high-density contrast could be generated 
by the excitation of slow-mode waves due to the non-linear 
steepening of a fast-mode wave. In LFH05, we explored how 
these effects were modified by the dissipation due to ion- 
neutral friction. The present paper is an extension of the 
work in FH02 to two dimensions. We find that the mecha- 
nism described in the previous paper works even better in 
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two dimensions in the sense that the density contrasts are 
both larger and persist for longer. 

The governing equations and the initial condition are 
given in Sect. [5] We then discuss the results for a specific 
numerical model (Sect. [3J and examine how these depend 
on the model parameters (Sect.|lJ. Finally, we conclude the 
paper by relating our results to observations of clumps and 
dense cores in GMCs (Sect. 01. 



2 NUMERICAL CALCULATIONS 
2.1 Governing equations 

While FH02 considered only variations in the :r-direction, we 
also allow changes in the y-direction. Since we are only inter- 
ested in magnetosonic waves, there is no need to introduce a 
z-component of either the magnetic field or the velocity. The 
MHD equations for an isothermal gas can then be written 
in the form 
dU | dF | dG _ Q 
dt dx dy 
where 

U = [p,pV X ,pVy,B X ,By] t 

is a vector of conserved variables, and 



F = [pv x ,pvl+p g + B 2 /2-Bl,pv x v y 

G = [pVy, pVyV X - ByB X , pVy + Pg + B 2 /2 



B X By,0, V X By 
B 2 y ,VyB x 



VyB x \ l 
VxBy^f 



are the corresponding fluxes. Here we assume that p g — 
pa 2 , where a is a constant sound speed. These equations 
w ere solved using the second order upw ind scheme described 
in lFalle. Komissarov fc .loarderl l)l998fl . excep t that we use 
the di vergence cleaning technique described bv lDedner et alJ 
(2002). 
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2.2 Initial state 

As in FH02, we are interested in the generation of density 
inhomogeneities by fast-mode waves propagating in the pos- 
itive ^-direction. Such a wave propagates with speed 

Vx + Cf, 

where the fast magnetosonic speed, Cf, is given by 



B 2 /p + a 2 + y/(B 2 /p + , 



4a a m/p 



C '=2 

The form of such a wave is most conveniently written 
in terms of the primitive variables 



P = [p,v x ,v y ,B x 



We then have 
9P 



dx 



OC Vf 



where 



P, c f, 



CfBxBy 

A, 



0. 



A/ 



with A 



2 

pc f 



Bm 



is the right eigenvector for a fast-mode wave propagating in 
the positive x direction for the system written in terms of 
the primitive variables, P. Note that we assume that the 
fast-mode wave is non-degenerate i.e. neither B y nor B x 
vanish. 

Like FH02, we use the above expression to superpose 
a small amplitude sinusoidal fast-mode wave propagating in 
the positive s-direction on a uniform background at rest. In 
order to make the flow two dimensional, we also introduce 
a shift with respect to the x-axis that depends on y. The 
initial state is given by 

' 2tt 



P(x,y,Q) = P + Ai sin 



r 2tt 
I IT 



X + j42Sin (l?) 



r/, (1) 



where A\ and L\ are the amplitude and wavelength of the 
fast-mode wave and A2 and L2 the amplitude and wave- 
length of the phase shift. The homogeneous background 
state Po is as in FH02: 

p=l, p = T c p, v = 0, B m = l, B y = aB x . 



Here T e (= 



is the equilibrium temperature of the 



isothermal gas. As in FH02, we use T c = 10 (or a = 
0.0316). Since the amplitude is small, the fast-mode 17 is 
calculated using the values of the unperturbed state and is 
normalised so that the amplitude of v x is A\. 

Density inhomogeneities are generated for all values of 
a. However, the highest contrasts are obtained when the ini- 
tial wave propagates at a modest angle to the field (FH02). 
In our calculations we therefore limit ourselves to a = 0.25. 
We set Ai = 0.05 since the transverse magnetic field changes 
sign for larger amplitudes . The non-linear steepening of the 
wave then leads to t he formation of an interme diate shock, 
which is unphysical iFalle fc Komissarovll200lT) . We adopt 
A2 = 1 for the amplitude of the phase shift. Other values 
for the model parameters are being considered in Sect. 2] 

Note that in the initial state, P(x,y,0), B x is constant, 
while By is a function of y. This means that the divergence 
constraint V.B = is violated. However, for the adopted pa- 
rameters (see also below), we find that V.B < 10~ 5 , which 
is of the same order of magnitude as the divergence errors 




Figure 1. The density at t = 2 000. The contour lines are added 
to show that two regions can be distinguished. The central region 
lies approximately within 250 < y < 750. 



that arise in numerical MHD calculations due to discretisa- 
tion of the e quations. As we inclu de the divergence cleaning 
algorithm of iDedner et alJ ll2002F) to stabilize the numerical 
scheme, the initial violation of the divergence constraint is 
damped out. 

Our initial state has an advantage over other models 
which introduce the two-dimensionality in a more natural 
way, i.e. a simple y-dependence of the primitive variables 
that allows an unambiguous interpretation of the initial 
stages of the simulations. A small density variation or a 
magnetic field B x that varies slightly with y introduces a 
complex y-dependence for the primitive variables (through 
Cf and Af). Calculations done with these alternative ini- 
tial states, however, produce similar results to those of our 
model and hence provide a posteriori justification for our 
choice. 

The computational domain is ^ x L\ and ^ 
y L2 (Li — L2 = L — 1000) with periodic boundary 
conditions. We use a uniform grid with 1 000 x 1 000 mesh 
points. 



3 RESULTS 

Initially the solution is almost one dimensional: the gra- 
dients in the {/-direction are of the order L/2ttA2{~ 100) 
smaller than the gradients in the x-direction. Although these 
variations are small, they cannot be neglected. From Eq. Q 
we see that the gradients in the y-direction change sign when 
y — L/4 = 250 and y — 3L/4 = 750 (for any given value of 
a;). This means that the domain can be divided into two re- 
gions with different temporal evolution: a central part with 
250 < y < 750 and the rest of the domain. In Fig. Q we 
can clearly see these two regions. Note that at this stage the 
variations in the y direction are of the order of a few percent. 
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Figure 2. (a) The growth of the density perturbation at y = 250. 
The times are t = 3 X 10 3 , 2 - t = 3.5 X 10 3 and 3 - t = 3.6 X 10 3 . 
(b) Decay of the density perturbation at y = 250. The times are 
1 - t = 3.7 x 10 3 2 - t = 4 x 10 3 and 3 - t = 4.5 x 10 3 . 



Figure 3. The maximum density measured in the computational 
domain for both ID (1; from FH02) and 2D (2) as a function of 
time. 




After 3 periods (i.e. t ~ 3 000) the fast-mode steepens 
into a shock, thereby exciting slow-mode components. FH02 
showed that these slow modes produce inhomogeneities with 
large density contrasts. Since there are as yet no significant 
variations in the y-direction, these high-density perturba- 
tions aligned with the j/-axis. Figure shows the growth of 
one such high-density region and its subsequent decay. Up to 
t ~ 5 000, the maximum densities in the computational do- 
main do not differ significantly from those in ID (see Fig.|3J. 
The density profiles of the perturbations are also similar to 
the ID ones (see Fig. Hand FH02). 

We can also see from Fig. that the density pertur- 



Figure 4. The density at t = 5 500. The density perturbations are 
bounded by curved shocks. The density is highest in the central 
region where the shocks are moving slowly. 



bation is bounded by shocks as it expands. In contrast to 
the ID case, the speed in the laboratory frame, v s , of these 
shocks is not constant along the shock front. The upstream 
fluid speed in the frame of the shock, v u , is roughly constant 
and equal to A\c a (FH02), where c a — \B X \/ p 1 ^ 2 is the ini- 
tial Alfven speed, whereas the upstream speed in the rest 
frame, v p , varies considerably along the shock front. Since 
v s — v p — v u , it follows that v s must also vary along the 
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shock front i.e. the shock front is curved. This effect can be 
seen clearly in Fig. [I] 

An important consequence of the curved shock front 
is the formation of clumps, rather than dense sheets. The 
dense regions are bounded by slow shocks that propagate 
predominantly in the x-direction. Regions that are bounded 
by slowly moving shocks therefore retain a high density for 
longer than those bounded by faster shocks. At later times 
this produces variations in the y-direction that are larger 
than those that occurred before the initial decay. However, 
there is little motion in the {/-direction because this is pre- 
vented by the magnetic field. Once a clump forms, its long 
axes is orientated perpendicular to the magnetic field lines, 
but its aspect ratio is not particularly large (see Fig. |SJ . 

Although the average density within clumps is nearly 
constant, the maximum density can have large fluctuations. 
As seen from Fig. the maximum densities are a factor 
of two higher than obtained in ID calculations. The initial 
fast-mode is still present and interacts with the clumps and 
its bounding shocks, thereby generating slow-mode waves. 
Figure |S]shows a measure of the amplitude of fast and slow- 
mode waves as defined in FH02. It shows clearly the cor- 
relation of the slow-mode waves with the density inhomo- 
geneities. The slow modes are thus responsible for the frag- 
mentation of the clump into high-density contrasts. Figure||j] 
shows a close-up of a clump with dense substructures. Such 
fragmentation is, of course, not possible in ID. 

The dense substructures in the clump are thus tran- 
sient and have rather short lifetimes. They are generated 
because slow-mode waves are continuously excited by the 
fast-mode wave interacting with the clump. This explains 
why the time-scales of the density fluctuations are of the 
order of the time for a fast-mode wave to cross the compu- 
tational domain, i.e. t = 1 000. As this wave dissipates, the 
amplitude of the slow modes decreases, resulting in lower 
density contrasts (see Fig. |5J. Furthermore, the clump it- 
self disperses. The 2D results show that a clump persists for 
a period of time that exceeds the lifetimes of the density 
perturbations in ID by a factor of 10. 

The ratio, /3, of the gas pressure to magnetic pressure in- 
creases in our simulations from its initial value of (5 m 0.002 
to P f» 0.2. While the gas pressure can become a factor 
of 100 larger than initially since it scales like the density, 
the magnetic field does not vary much in strength. When ft 
approaches unity, the process of producing high-density con- 
trasts by slow modes ceases to be effective. There is there- 
fore an upper limit to the maximum density that can be 
produced. 



4 PARAMETER STUDY 

In the previous section, we focused on one specific numer- 
ical case with a set of specific parameters. As the results 
are parameter-dependent, we briefly discuss the influence of 
each parameter. 




0.0 0.2 0.4 0.6 





a : The fast-mode wave is most efficient in exciting slow- 
mode components for modest values of the angle, a, between 
its direction of propagation and the magnetic field (FH02). 
When a is too large, the shock is not effective in produc- 
ing slow-mode components, whereas there is little steepen- 



Figure 5. (a) Similar to Fig. 131 but now for a uniform grid with 
100 X iOO mesh points up to t = 100 000. (b) The maximum 
value for the amplitude-measure of the slow-mode wave in the 
computational domain (for definition see FH02). (c) Similar to 
(b), but for the fast- mode wave. 
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Figure 6. Hierarchical structure in the density at t = 12 000: 
a number of dense cores reside in a larger clump. The magnetic 
field is represented by the solid lines. 

ing when a is too small. Moreover, the fast-mode wave is 
also responsible for the excitation of slow modes within the 
clumps. Since it is the slow modes that generate the den- 
sity inhomogeneities, this means that models with a 7^ 0.25 
produce lower density contrasts. 

Ai and a : Both parameters affect the maximum den- 
sity produced by the non-linear steepening of the fast-mode 
wave. FH02 showed that the maximum density is roughly 
proportional to A\/a 2 . 

While the time-scale for the formation of the substruc- 
ture is independent of a, it depends on Ai in two ways (see 
FH02). The wave amplitude determines the time for the fast- 
mode wave to steepen into a shock and also the time-scale 
for the subsequent decay. A lower value of Ai means that 
dense cores arise later in the simulations. 

A2 : The phase-shift of the fast-mode wave is responsible 
for the variations along the y-axis. The amplitude of this 
shift determines the time at which the simulations cease to 
be ID. Figure Q shows the time evolution of the maximum 
density for an amplitude of A2 = 0.1. As the initial vari- 
ations are a factor of 10 smaller than for A2 = 1, it takes 
longer to develop any significant deviation from ID. This 
also means that the formation of substructure is deferred. 
Note that the maximum density that can be attained is in- 
dependent of A2 (see Fig. |3] and Fig. 0. 

Li and L2 : The wavelengths of the fast-mode wave and 
the phase shift determine the dimensions of the clump. Its 
extent in the ^-direction is roughly comparable with the 
size of density perturbations produced by excitation of slow 
modes. This means that it is ~ Lipo/Pmax, where po is the 
initial density and p max the maximum density. The size in 
the y-direction is set by the wavelength of the phase shift. 
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Figure 7. Similar to Fig.|3 but now for A2 = 0.1. 

As mentioned in Sect. [3] a sinusoidal phase shift divides 
the domain into two separate regions so that the y-extent is 
roughly L2/2. 



5 DISCUSSIONS AND CONCLUSIONS 

We have extended the FH02 model which describes the gen- 
eration of density inhomogeneities by MHD waves in a low-/? 
environment, to two dimensions. We have followed the evo- 
lution of a fast-mode wave propagating in the positive x- 
direction with a small, but non-negligible, variation in the 
y-direction produced by a y-dependent shift along the a;-axis. 
Like FH02, we find that the formation of density inhomo- 
geneities is associated with the excitation of slow modes by 
the fast-mode wave. 

An important result is that the process of generating 
structure by slow modes works on different length-scales. 
While the larger structures arise due to the non-linear steep- 
ening of the fast mode, there are smaller structures gener- 
ated by the interaction between the fast-mode wave and 
the larger structures. A limiting factor, however, is that 
P increases within each level of substructure. When /3 ap- 
proaches unity, the mechanism ceases to be effective and 
produces no further substructur es. j3 is probably about unity 
in dense cores 

Although we made no attempt at reproducing any ob- 
servational data, we can investigate the relevance of our 
simulations to observations of clumps and dense cores in 
GMCs. The structure of molecular clo uds is often de- 
scrib ed as filamentary or clumpy (e.g. iBlitz fc Williamsl 
119971) . For example, observations of the p Oph A 
core show an arc-shaped ridge with 4 embedded cores 
jAndre. Ward- Thompson fc Barsonvl 119931. . Similar struc- 
tures are obtained in our simulations. However, the shape of 
the clump depends on the properties of the initial fast-mode 
wave, i.e. its wavelength, and on the length-scale of per- 
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pendicular variations (in our calculations determined by the 
sinusoidal phase-shift of the fast-mode wave). Other shapes 
of clumps can then be reproduced by MHD waves with dif- 
fere nt properties. 

iBlitz fc Stark! (1986) found that the clumps in the 
Rosette Molecular Cloud, an archetypal GMC, are an order 
of magnitude denser than the mean cloud density. Dense 
cores within the clumps typically have a density that ex- 
ceeds the mean density by two orders of magnitude, but 
peak densities in pre-stellar cloud cores can range to values 
that exceed t he mean density by 4 orders of magnitude (e.g. 
iLarsonlEoO^I . In our simulations, we find that the average 
density contrasts in the larger structures are in excess of 10 
times the average background density of 1. For the substruc- 
tures, i.e. individual dense cores, this increases to 100. We 
do not find the observed peak values, but this can be read- 
ily explained by the fact that we have restricted ourselves 
to modest amplitudes of the initial fast-mode wave. Also 
we have not included self-gravity in our model. The inclu- 
sion of this effect will be investigated in a subsequent paper. 
However, we expect self-gravity only to act as an amplifier 
and not to have an influence on the fragmentation of the 
molecular clo ud. 

Recently, iGarrod et alj i2005l) investigated the chem- 
istry of transient dense cores within molecular clouds. Using 
the results of FH02 to describe the growth and decay of the 
dense cores, they find a significant enhancement in chem- 
ical compositions of the cores as they disperse. This leads 
eventually to chemical enrichment of the molecular cloud. In 
our simulations, however, high-density clumps survive much 
longer than they do in FH02. Our results thus have seri- 
ous implications for the chemical evolution and richness of 
clumps and dense cores. 
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